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SUMMARY 


A survey of methods for sensitivity analysis of the algebraic eigenvalue problem 
for non-Hermitian matrices is presented. In addition, a modification of one method 
based on a better normalizing condition is proposed. Methods are classified as Di- 
rect or Adjoint and are evaluated for efficiency. Operation counts are presented in 
terms of matrix size, number of design variables and number of eigenvalues and 
eigenvectors of interest. The effect of the sparsity of the matrix and its derivatives is 
also considered, and typical solution times are given. General guidelines are estab- 
lished for the selection of the most efficient method. 


Introduction 

The behavior of many physical systems is completely determined by the 
eigenvalues of the system. Variations in parameters lead to changes in these 
eigenvalues and hence in response characteristics. Thus derivatives of eigenvalues 
and eigenvectors are of immense interest in several fields of physical sciences and 
engineering. 

The derivatives (or synonymously, sensitivities) are of interest for a variety of 
uses. Design optimization is intimately connected with sensitivity analysis and the 
cost of calculating derivatives is the dominant contributor to the total cost in many 
optimization procedures. Most optimization algorithms require many analyses of the 
system and derivatives can be effectively used to approximate the eigenvalues and 
eigenvectors of a modified system and thus reduce the cost of reanalysis, especially 
in large systems. In addition, derivatives are very useful in design trend studies and 
for gaining understanding of and insight into the behavior of physical systems. 
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Finally, derivatives of eigenvalues are valuable in calculating the statistics of 
eigenvalue locations in probabilistic analyses. 

The application of derivatives is not restricted to design-oriented activities. 
Sensitivity analysis is also playing an increasing role in determining the analytical 
model itself. In the areas of system identification and analytical model improvement 
using test results, sensitivity analysis is of growing importance. Much recent work 
in these fields is directly dependent on the calculation of eigensystem derivatives. 

It has been found in certain cases that second order derivatives are effective in 
improving accuracy of approximations[1-7] and efficiency of design[3,8,9]. 
Eigenvalues are usually non-linear functions of design parameters and a second or- 
der approximation offers a much wider range of applicability compared to the first 
order approximation. Intermediate variables which may improve the quality of first 
order approximations are not generally available for eigenvalue approximations. 
Also, some optimization algorithms require second order derivatives, and first order 
derivatives of optimal solutions require second order derivatives of constraints[7]. 
The use of second derivatives can also greatly reduce the number of reanalyses re- 
quired for the convergence of an optimization procedure^, 10]. Further, in certain 
optimization algorithms, second order approximations for eigenvalue constraints can 
drastically relax the move limits, thus achieving a nearly optimum trajectory, and can 
virtually eliminate the need for trial and error adjustment of move limits, thus im- 
proving the performance of the optimizer[10]. Looking at another aspect, in problems 
where instabilities are to be avoided, a first order calculation may completely fail to 
detect instabilities^]. References [11,12] also offer examples of the usefulness of 
second order derivatives. 

The problem of calculating the derivatives of symmetric and Hermitian 
eigenproblems is relatively simple and solution procedures are well-established, 
e.g. [13-17], However, many physical problems give rise to non-self-adjoint formu- 
lations and thus lead to general matrices. An important example is aeroelastic sta- 
bility which requires the solution of eigenproblems with complex, general and fully 
populated matrices. General systems are also obtained in damped structural sys- 
tems anc[ in network analysis and control system design where the eigenvalues are 
usually called poles. This study presents a comparative analysis of the various 
methods available for calculating the derivatives of the general eigenproblem and 
propose some modifications to existing techniques. A considerable amount of liter- 
ature is available on the subject and a comparative analysis of the various methods 
will be of value for selecting the most efficient technique for a particular application. 
The purpose of this paper is to summarize the more efficient techniques proposed 
so far and to establish guidelines for the selection of the appropriate method for a 
given problem. Only the essentials of these methods are presented with details re- 
ferred to the original references. Attention is restricted to the general eigenproblem 
and techniques that are useful only for the self-adjoint problem are not considered. 
The present discussion is limited to the case of distinct eigenvalues. 



Problem Definition 


The matrix eigenproblem is defined as follows: 

Au“> = X'V (1) 


and the corresponding adjoint problem is 

v<*> 7 A = \ {k V k)T (2) 

where A is a general complex matrix of order n and uW and vW are the k -th 
eigenvalue and right and left eigenvectors respectively. The superscript T denotes 
the transpose. All eigenvalues are assumed to be distinct. 


The matrix A and hence, uW and vW are functions of design parameter 
vector p with individual parameters denoted by Greek subscripts, e.g. p a . Deriva- 

<9A 

tives with respect to p a are denoted by the subscript , a e.g., — — = A a . All the 
design variables are assumed to be real. ^P a 


The well-known biorthogonality property of the eigenvectors is given by 
v ( ' )T u W = 0 iff / * j (3) 


Note that, the left hand side of eq. (3) is not an inner product as usually under- 
stood, since vW and/or u W may be complex vectors. The left eigenvectors of A are 
the right eigenvectors of A 1 and vice versa. 


Normalization of Eigenvectors 


The eigenvectors uM and vW are not completely defined by eqs. (1) and (2). A 
normalization condition has to be imposed to obtain unique eigenvectors. For brev- 
ity, let us consider only the normalization of the right eigenvector. A normalizing 
condition frequently imposed in the self-adjoint case is the following: 


U WV*> = i 


(4) 


However, it is not always possible to use eq. (4) for non-self-adjoint problems as 
u (k)T u W can e q Ua | zero 0 r a very small number causing numerical difficulties. This 
is true even if the matrix A is real. Unfortunately, considerable confusion exists in the 
literature regarding this point and several authors arbitrarily adopted eq.(4) as a 
normalizing condition for non-self-adjoint problems, e.g. [8, 9, 11, 18-21], In this re- 
spect, the formulations of these references are not rigorous for general matrices. 

One possible way to avoid the above difficulty is to replace eq.(4) by 
«,<*>«<*> = 1 (5) 


where superscript * denotes a conjugate-transpose. Eq. (5) is not prone to the diffi- 
culties of eq. (4) because uM*uM is always guaranteed to be non-zero. But, eq.(5) is 
not a complete normalizing condition as it does not render the eigenvector unique. 
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If u satisfies eq.(5), then w = ue ,c , where / == >/ — 1 and c is an arbitrary real num- 
ber, also satisfies eq.(5). 


Another normalization condition, inspired by the biorthogonality property of the 
left and right eigenvectors, is 

V wv*> = i (6) 


Eq.(6) also does not render the eigenvectors unique. It must be emphasized that 
if the eigenvector is not unique, neither is its derivative. 


The normalization condition 



= 1 


( 7 ) 


is very attractive because it renders the eigenvectors unique and at the same time, 
the index m can be chosen easily to avoid ill-conditioning. Apparently, only 
Nelson[22] used this normalizing condition in obtaining the derivatives of 
eigenvectors. The index, m, may be chosen such that 

|ciJU ) l = mpx | u} k) | (8) 


Another choice for m, used by Nelson[22], is 
lu^l |i$>| = max | uj k) | (9) 


The nature of uncertainty of the derivative of the eigenvector is of some interest. 
Without a normalizing condition, an eigenvector is uncertain to the extent of a non- 
zero constant multiplier. The derivative of an eigenvector is uncertain to the extent 
of an additive multiple of that eigenvector. To show this, let uM be an eigenvector 
so that wW = cuW is also an eigenvector. Then, if p a is a design parameter, 

d(c u W ) 

d Pa d Pa 

where d = (dc/dp a ) is arbitrary. In practice, the constant d depends on the way the 
eigenvectors uW and wW are normalized. 



<?!»<*> 

d Pa 


+ du 


W 


( 10 ) 


Methods of Calculation 

The various methods of calculating the derivatives of eigenvalues and 
eigenvectors can be divided into three categories: 

1. Adjoint Methods, which use both the right and the left eigenvectors. 

2. Direct Methods, which use only the right eigenvectors. 

3. Iterative Methods, which use an iterative algorithm that converges to the re- 
quired derivatives. 
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Adjoint Methods 


The first expressions for the derivatives of eigenvalues of a general matrix seem 
to have been derived by Lancaster[23J. Considering only a single parameter, 
Lancaster obtained the following expressions for the first and second derivatives of 
an eigenvalue: 


W = 


V W 7 A a u<*> 


vWy*) 


v WT A u 

Jk) = v M . cta u 

' a “ v W7 u W 


w 


n 

+ 21 
j= 1 
j*k 


(v^ 7 A a u^) (v 


W r A 0 uW) 


- xW)(v<*)y*)) (V wr u w ) 


( 11 ) 

( 12 ) 


An expression corresponding to eq. (11) for a generalized quadratic eigenvalue 
problem was obtained by Pedersen and Seyranian[24]. 

Morgan[25] obtained an expression for the derivative of an eigenvalue without 
requiring the eigenvectors explicitly. His expression is equivalent to 

. _ trace of {[adj(A - X w I)]a „} 

trace of adj(A - X ( ^l) 

The corresponding expression for derivatives with respect to matrix elements was 
derived by Nicholson[26], 


It can however be shown that[27] 
adj(A - = t k u w v^ T 

where t k is a constant and that[28] 

trace of {[adj(A - X W I)]A = t k v {k)T A a u w 
trace of adj(A - )J k h) = t k v^ T 


Thus, in the computation of adj^A — both right and left eigenvectors are 

implicitly computed, in view of eq. (14). Eqs. (15) also show that Morgan's eq. (13) is 
equivalent to Lancaster's eq. (11). Woodcock[29] also obtained formulas involving 
the adjoint matrix for the first and second derivatives of eigenvalues. An operation 
count shows that calculation of the adjoint matrix is several times more expensive 
than the explicit calculation of right and left eigenvectors so that Lancaster's formula 
is preferable to formulas requiring the adjoint matrix. This conclusion is also sup- 
ported by sample computations[30j. In addition, although eq. (13) was used satis- 
factorily for small problems[31,32], numerical difficulties were reported for 
reasonably large problems[33]. Woodcock's formula for the second derivative of the 
eigenvalue requires a partial derivative of the adjoint matrix and this is so compli- 
cated that Woodcock himself recommends the finite difference method. Formulas 
due to Morgan and Woodcock are not therefore considered in the following. 
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To obtain the second derivatives of eigenvalues, the first derivatives of left and 
right eigenvectors are calculated either implicitly[9,11,23] or explicitly[1,8,12,34,35]. 
Since the eigenvalues are assumed to be distinct, the first derivatives of eigenvectors 
can be expressed as 


u !a = Zckja " 01 and 


,W = v 


i dn,y 


Rogers[36] obtained the coefficients c k j a and d k j a as 


v^A a u<*> 


rn 


d kja = 


(x w - xW) 


( i)\ v U)T u U ) 


vW r A a u» 

.(*) - 1 0‘)\,.(/)T 11 (/) 




V w U 


k +J 


k *j 


It can be observed that 

d k ;„ = - C;, ra - (19) 

kja Jka V (/)T U 0) 

Reddy[37] derived an equivalent expression for the response derivative by casting 
the derivative as the solution of a forced response problem for the same system. 

Note that, in view of eq. (10), the coefficients c kka and d kka in eq. (16) are arbitrary 
and depend on the normalization of the eigenvectors. For example, if eq. (7) is used 
to normalize the right eigenvectors, then 

c kka = ~ ^ c kja u m (20) 

7=1 

i*k 

and if eq. (6) is used to normalized the left eigenvectors, then 


d kka ~ c kka 


It has been proposed[38,39] that the eigenvector derivative be approximated by 
using less than the full set of eigenvectors in the expansion of eq. (16) so that the 
evaluation of eigenvector derivative by Adjoint Method could become cheaper. This 
variant of Adjoint Method has received mixed reports in the literature[22,38]. The 
quality of such an approximation is difficult to assess beforehand and the selection 
of eigenvectors to be retained in the expansion is problem dependent. It has not 
been considered in this work because a meaningful comparison with other methods 
cannot easily be made. However, this consideration should not be ignored while 
implementing the sensitivity calculations for particular problems. 
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The expressions for the second derivatives of eigenvalues were obtained by 
Plaut and Huseyin[35]. For the sake of simplicity in expressions, let us assume, 
without loss of generality, that the left eigenvectors are normalized as in eq. (6). 


Plaut and Huseyin[35] obtained the second derivatives of eigenvalues with respect 
to uncorrelated parameters p a and pp as 


^ - v« r A ap u<‘> + v<> a u^ + v^A a u<*> 


( 22 ) 


which can be equivalently written, without involving the derivative of the left 
eigenvector, as 


*■% = v WT A op u W + v Wr (A a - \ a l)u^ + v WT (A p - X pDu^ 
For a diagonal second derivative, eq.(22) is simplified to 

= v <*>r A u (« + 2v (*)r A „(*) 

, aa v , aa u M , a u , a 

Eq. (22) can be rewritten using eqs. (17) and (18) as 

= v<* |r A + .£ + c m d k]a ) 

)*k 


(23) 


(24) 


(25) 


Crossley and Porter[1,40] derived similar expressions for derivatives with re- 
spect to the elements of the matrix. Expressions for AMh order diagonal derivatives 
were derived by Elrazaz and Sinha[5]. 

In calculating the derivatives using eqs. (11), (16)-(25), 

• the first derivative of an eigenvalue requires the corresponding right and left 
eigenvectors. 

• the first derivative of an eigenvector requires all the left and right eigenvectors. 

• the second derivative of an eigenvalue requires the corresponding right and left 
eigenvectors and their first derivatives. 


Direct Methods 

The second category comprises methods that evaluate the derivatives using only 
the right eigenproblem. Direct Methods typically involve either the evaluation of the 
characteristic polynomial or the solution of a system of linear simultaneous equations 
without requiring all the left and right eigenvectors. Methods requiring the evalu- 
ation of the characteristic polynomial and the derivative of the determinant[33,41] are 
0(n 5 ) processes while other methods considered here are at most 0(n 3 ) processes. 
In addition, the determination of the characteristic polynomial is, in general, an un- 
satisfactory process with respect to numerical stability, even when all the 
eigenvalues are well-conditioned[42]. While numerically stable algorithms have been 
proposed for evaluation of the characteristic polynomial[43], the computational ex- 
pense still seems to be formidable. Hence, we do not consider these methods. 
Methods requiring the solution of a system of equations have the particularly attrac- 
tive feature that the coefficient matrix needs to be factored only once for each 
eigenvalue regardless of the number of parameters and the order of the derivatives 
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required. Thus, they are very useful in applications where higher order derivatives 
are required. 

The earliest method in this class is due to Garg[18] who obtained the first deriv- 
atives of the eigenvalue and the eigenvector by solving two systems of (n + 1) 
equations each in the real domain, without requiring any left eigenvectors. However, 
his formulation involves several matrix multiplications. Rudisill[19] proposed a 
scheme in which only the corresponding left and right eigenvectors are required to 
calculate the first derivative of the eigenvalue and the eigenvector. This was refined 
by Rudisill and Chu[20] to avoid calculating the left eigenvectors altogether. Solution 
of a system of only ( n + 1) equations is required (though in the complex domain) to 
obtain the first derivatives of eigenvalue as well as eigenvector. Extension to higher 
order derivatives is straightforward. Cardani and Mantegazza[21] proposed solution 
methods of the same formulation for sparse matrices and extended it to the quadratic 
eigenproblem. 

One weakness that is common to all the above formulations that do not require 
left eigenvectors[18-21] is that they rely on the normalization condition given by eq. 
(4), which is unreliable for general eigenproblems as discussed earlier. 

Nelson[22] circumvented this difficulty by using the normalizing conditions 
v ( * )r uW-1 and „W -1 


However, the formulation of Rudisill and Chu is superior to Nelson's formulation 
in that it does not require any left eigenvectors. 

In this paper, we propose a variation of the Rudisill and Chu formulation which 
does not rely on the questionable normalizing condition of eq. (4) and at the same 
time requires no left eigenvectors. 


Differentiating eq. (1), we get 

Au^A/) = xV>+Jl!i ) u W 

which can be rewritten in partitioned matrix form as 

[A -JL«I | -A. 0 U<‘> 


(26) 


(27) 


Now, we impose the normalizing condition of eq. (7). Differentiation of eq. (7) 
yields, 


= 0 

m, a 


(28) 


Because of eq. (28), the m-th column of the coefficient matrix in eq. (27) can be 
deleted. Eq. (28) also reduces the number of unknowns by one so that eq. (27) is now 
a system of n equations in n unknowns. Eq. (27) is rewritten as 

B yi = r (29) 
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where 


B = 

r a _ i Mi | _ ..Mi 

t M K 1 1 u 1 m-th column deleted 



f<\ 


Vi = 

\ X {k) I 

>■ , a J with m-th element deleted 


r = 

1 

P 

5 

(30) 


To get second derivatives, differentiate (27) with respect to pp and get, 
(A - - u“\ ap - - A, ap u<‘> - (A. - X;*>)u^ p > 

- (A, „ - 

or, in partitioned matrix form, 

r J A 

[A - | - - A, apU «‘» - ,A „ - 

- (A, p -*><*> 

Following the same reasoning as before, eq. (32) is written as 
By 2 = s 

where 


V2 = 




(31) 


(32) 


(33) 


J with / 77 -th element deleted 


s = 


- A „u W - 
, aP u 


(A 


- iWiJ*) - 


- (A „ - ^'l)u»» 


(34) 


Note that, if is a simple eigenvalue of A and if u$ #= 0, then the matrix A is 
of rank (n — 1) and the m-th column that is deleted is linearly dependent on the other 
columns. Hence the matrix B is non-singular. The matrix B will also be well- 
conditioned if uffl is the largest component in the eigenvector uW and the matrix A 
is itself not ill-conditioned. The vectors y-| and y 2 can be obtained by standard sol- 
ution methods. If the matrix A is banded or if the derivatives of both right and left 
eigenvectors are required, it may be more efficient to use a partitioning scheme as 
described in the appendix. 


In summary, we note that, in calculating derivatives by Direct Method, 

• left eigenvectors are not used. 

• a complete solution of the eigenvalue problem is not required, if the derivatives 
of only a few of the eigenvalues and eigenvectors are sought. This is in contrast 
to the Adjoint Method which requires all the left and right eigenvectors to cal- 
culate the first derivative of any eigenvector. 
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• calculation of any derivative requires the solution of a system of linear 
equations. 

• only one matrix factorization needs to be performed for all orders of derivatives 
of an eigenvalue and its corresponding right and left eigenvectors. 


Iterative Methods 

Andrew[44] proposed an iterative algorithm to calculate the first derivatives of 
eigenvalues and eigenvectors. This algorithm is a refined and generalized version 
of the iterative scheme developed by Rudisill and Chu[20]. Except for the dominant 
eigenvalue, the convergence of this algorithm seems to be very much dependent on 
the choice of the initial values for the derivatives. To be efficient for non-Hermitian 
matrices, this iterative method requires a complex eigenvalue shifting strategy which 
is not easy to implement. Hence this method is not considered. 


Relative Computational Cost 

In this section, we compare the efficiency of calculating the derivatives of 
eigensystems as a function of the size of the matrix n , number of design parameters 
m and number of eigenvalues of interest /. 

To start with, let us consider the operation counts (flops) for the Adjoint Methods 
given by eqs. (11),(16)-(25) and the Direct Methods given by eqs.(29)-(34). They are 
summarized in Table 1. It should be noted that the operation counts represent an 
estimate of the actual number of operations performed by a solution routine and in- 
clude only the most significant terms. The actual number of operations will vary 
slightly depending on programming details. The effect of the sparsity of the matrix 
derivative A a is modeled by the parameter k, defined such that the the number of 
operations in evaluating the product A a u is equal to K/? 2 (that is, k = 1 corresponds 
to a full A H ). 

The eigenvalues are calculated using the EISPACK subroutine package [45] by 
first reducing the matrix to upper Hessenberg form using unitary similarity transf- 
ormations and then applying the QR algorithm. The number of operations and CPU 
time for calculating the eigenvalues are not counted in the following results as they 
are not relevant in comparing the methods to calculate the derivatives. 

The right eigenvectors are calculated by inverse iteration on the same upper 
Hessenberg matrix used for calculating the eigenvalues and back transformation us- 
ing standard subroutines in the package EISPACK. The corresponding operation 
count is given in Table 1. For the calculation of left eigenvectors, it is important to 
note that there is no need to repeat the process with the transposed matrix. The left 
eigenvectors are obtained cheaply using forward substitution in place of backward 
substitution in the inverse iteration process. There is also no need to repeat the 
matrix factorization. A subroutine was written to calculate the right and left 
eigenvectors in this manner and the corresponding operation count is given in Table 
1 . 
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Table 1 gives the operation count of evaluating the individual steps. To obtain 
the number of operations involved in evaluating the derivatives, we must add the 
operation counts for all the steps required in the calculations. These counts are 
given in the following discussion. 


CPU Time Statistics 

In the following tables, computational cost for the calculation of the first and 
second derivatives of eigensystems are compared for matrices of order 20, 40 and 
60. The CPU time statistics are obtained on the IBM 3084 computer using the 
VS-FORTRAN compiler with no compiler optimization. The ratio of operation 
count(OC) and CPU time for various operations is about 10 5 operations per CPU 
second with a variablity of 27 percent. 

The matrices are generated for the dynamic stability analysis of a compressor 
stage rotor with mistuned blades. The geometric and structural parameters of the 
rotor and formulation and method of analysis are the same as those of NASA Test 
Rotor 12 described in reference[46] except that the number of blades and the 
torsional frequencies are varied. The torsional frequency values are selected ran- 
domly from a population of mean 1.0 and standard deviation 0.01. The standard de- 
viations of the actual samples are slightly different. 


Calculation of First Derivatives of Eigenvalues Only 


Operation Count 


Adjoint Method 
Direct Method 


/(-^-n 2 + Kmn 2 ) 
/[-j- + («c + 1)mn 2 ] 


It is clear from the operation count that the Adjoint Method, which is an 0(n 2 ) 
process, is superior to Direct Method, an 0(n 3 ) process, for large n. The number of 
design variables and the number of eigenvalues of interest have no bearing on this 
conclusion. As the order of the matrix increases, the Direct Method becomes more 
expensive. For example, for 5 design variables and 10 eigenvalues of interest, the 
CPU time for the Direct Method is 2.3 times more expensive than for the Adjoint 
Method for n = 20, and for n = 60, the ratio is 3.0. 


187 


Calculation of First Derivatives of Eigenvalues and Eigenvectors 


Operation Count 


Adjoint Method 


Direct Method 


-^-n 3 + lmn 2 { k + 2) 
— f lmn 2 ( k + 1) 


When the derivatives of both eigenvalues and right eigenvectors are required, 
the choice of method is dependent on the values of / and m. When very few 
eigenvalues are of interest, the Direct Method is cheaper. When many eigenvalues 
are of interest, the Direct Method is more expensive than the Adjoint Method. How- 
ever, this effect of the number of eigenvalues of interest is less significant when the 
number of design variables is large. As the number of design variables increases, 
the Direct Method becomes more competitive, even when all eigenvalues are of in- 
terest. For a 60 x 60 full (k = 1) matrix, this is illustrated in Figure 1. 

The operation count shows that the computation by Adjoint Method of 
eigenvector derivative, which is necessary for the second derivative of eigenvalue, 
is an O (n 3 ) process and is more expensive than the computation of the eigenvector 
itself which is an 0(n 2 ) process using the procedure described above. This fact is 
significant as some authors have stated the opposite[2,3]. 
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Calculation of First and Second Derivatives of Eigenvalues only 


Operation Count 

Adjoint Method 
Direct Method 
Direct-Adjoint Method 


-~/7 3 + (k + 1)mn 3 + / n 2 K 

lJ T + l ( 2 ) + 1 ) 

n 2 K + lmn 2 (2K + 1) 


The Direct-Adjoint Method denotes the calculation of the eigenvector derivatives 
by the Direct Method and the eigenvalue derivatives by the Adjoint Method. The 
third term in the operation count for the Direct-Adjoint Method is significant only 
when m is small. From the operation count, it is seen that the Direct-Adjoint Method 
is always cheaper than the Direct Method. Hence, the choice is between the Adjoint 
Method and the Direct-Adjoint Method. Here, considerations similar to those of the 
last section hold and the choice of method depends on the values of / and m. When 
few eigenvalues are of interest, the Direct-Adjoint Method is cheaper. When many 
eigenvalues are of interest, the Adjoint Method is superior. But this advantage of 
Adjoint Method diminishes as the number of design variables increases. This is 
again illustrated for a 60 x 60 full matrix (k = 1) in Figure 2. 


Concluding Remarks 

The normalization of the eigenvector needs to be properly related to its deriva- 
tive. In practice, this means that the derivative of the eigenvector is to be normalized 
before it is used, to conform to the normalization of the eigenvector itself. When the 
eigenvector is not normalized in a unique manner, its derivative cannot be evaluated. 
Fixing one of the components of the eigenvector is the best normalizing condition for 
computation of the derivative. The methods found in the literature are extended to 
apply to eigenvectors normalized in this manner. 

Various methods for calculation of derivatives of eigenvalues and eigenvectors 
are surveyed and classified as Direct or Adjoint. Adjoint Methods use both the left 
and the right eigenvectors whereas the Direct Methods use only the right 
eigenvectors. Their relative efficiency is evaluated as a function of matrix size, 
number of eigenvalues of interest and the number of design parameters. General 
recommendations are made for the cases when (a) eigenvalue first derivatives are 
required, (b) eigenvalue and eigenvector first derivatives are required, and (c) 
eigenvalue second derivatives are required. 
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When only eigenvalue first derivatives are required, the calculation of left 
eigenvectors is worth the expense as the Adjoint Method is shown to be superior to 
the Direct Method. When first derivatives of eigenvectors are also required, the de- 
cision is dependent on the problem size, the number of design variables and the 
number of eigenvalues of interest. When the first and second derivatives of 
eigenvalues are required, similar considerations hold. It is also shown that once the 
first derivatives of eigenvectors are calculated, the second derivatives of eigenvalues 
are calculated more efficiently by the Adjoint Method than by the Direct Method. 
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Appendix 


Modification of Direct Method for Banded Matrices 

Equations (29) and (33) can be written as 
(A X l)m-th column deleted 0 , a m-th row deleted — \ a a — (^1) 


Let uW be normalized so that u$ = 1 

Eq. (A1) is a system of n equations. Writing the m-th equation separately, we 
have, if the superscript (k) is omitted for notational convenience, 

Cx , a - \ a x = t (A2) 


and 

a m x , a — \ a — r m 

where 

^ — ^ ~ AJ)/n-th row and column deleted 
x , a — u , a m-th row deleted 
x = u m-th row deleted 
* ~ r m-th row deleted 

a m = m- th row of A with the m-th column deleted 
From (A3), 

\ a = a m x , a ~ r m 
From (A2), 


(A3) 


m 
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X >a = c V a X + t) 
Eliminating x a , we have 


where 

b m - [C r r'a m 

Proceeding in a similar manner for the left eigenvector, 

y,„ = (C T 1~ V „y + t,) 

where 


(A5) 


m 


(A 7) 


y, a v , a m - th row deleted 
y = v m- th row deleted 
t/ - fr/) m - th row deleted 
r/ being the appropriate right hand side. 

Thus the following procedure can be used to obtain the derivatives X a and u a . 

1. Form a LU decomposition of the matrix C. 

2. Solve b m = [C r ] - 1 a m by forward substitution. 

3. Calculate X a from (A6). 

4. Calculate x a from (A5) by backward substitution. 

5. Expand x a to u a setting u m a = 0. 

If the derivatives v a of the left eigenvectors are also required, only three further 
steps are needed. 

6. Calculate y a from (A7) by forward substitution. 

7. Expand y a to v a setting v m a = 0. 

8. Normalize v a appropriately depending on the normalization of v. For example, 
to obtain the derivative of the left eigenvector that satisfies the normalization 
condition of eq. (6), subtract (v r u a + v^ a u)v . 

The matrix C needs to be factored only once. Also, the matrix C retains the 
bandedness characteristics of the original matrix A. Furthermore, higher derivatives 
can be obtained by merely substituting an appropriate right hand side vector, r. 
However, higher order derivatives can suffer in accuracy because of accumulated 
round-off error. 
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The conditioning of matrix C needs some comment. Note that C is obtained from 
the singular matrix (A - A,M|) by deleting both the row and column corresponding to 
index m. Hence, for matrix C to be non-singular, one must make sure that the m-th 
row is linearly dependent on the other rows as well as that the m-th column is line- 
arly dependent on the other columns. In other words, C is non-singular iff 
u$ * 0 and vffl * 0 . If vffl is very small compared to the largest element in vM, 
steps 2 and 4 in the above procedure will give inaccurate results even if u$ is the 
largest element in uW. | n general, it is not possible to make a good choice for m 
without the knowledge of the left eigenvector. Since the calculation of left 
eigenvector by forward substitution is cheap, it is suggested that the left eigenvector 
be calculated and the index m be chosen as in eq.(9). This is the same criterion used 
by Nelson[22] and will assure as well-conditioned a matrix C as possible. 
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Table 1. Operation Counts 


Eigenvectors 


Operation 

Evaluation of right eigenvectors 
Evaluation of left eigenvectors 


Operation Count 
/( 2n 2 ) 

/( Y" 2 ) 


Adjoint Methods 


Operation Operation Count 

Evaluation of eq. (11) lmn 2 K 

Evaluation of eqs . (16),(17),(18) lmn 2 { k + 2) 

Evaluation of eq. (25) / n 2 k 


Direct Methods 


Operation Operation Count 


LU decomposition of matrix B 
Formulation and solution of eq.(29) 

Formulation and solution of eq.(33) 


l{n 3 /3) 

lmn 2 (K + 1) 


'( ™)n 2 (3K + 1) 


CPU TIME 
(SECONDS) 
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Figure 1. CPU Times for calculation of first derivatives of eigenvalues and eigenvectors for a 60 x 60 
matrix 
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Figure 2. CPU Times for calculation of second derivatives of eigenvalues for a 60 x 60 matrix 


